#include <iostream>
#include <fstream>
#include <chrono>
#include <ratio>
#include <iomanip>
#include <random>
#include <thread>
#include <functional>
#include <vector>
#include <string>
#include <CL/sycl.hpp>

using namespace cl::sycl;
using gauss_func = std::function<void(buffer<float, 2>&, queue&)>;

std::string data_path = "../data.txt";
std::string outputResult_path = "../outputResult_GEC_P.txt";

void matrixPrint(buffer<float, 2>& buf)
{
	host_accessor m{ buf ,read_only };
	auto range = m.get_range();
	for (int i = 0; i < range[0]; i++) {
		for (int j = 0; j < range[1]; j++) {
			std::cout << std::setw(16) << m[i][j];
		}
		std::cout << std::endl;
	}
}

// 先把数据读取到内存
void dataInit(std::vector<int>& Scale, std::vector<std::vector<float>>& Data)
{
	std::ifstream i_file;
	i_file.open(data_path, std::ios::in); // 读取数据
	if (!i_file.is_open())
	{
		return;
	}
	while (!i_file.eof()) // 文件成功打开
	{
		std::string data;
		while(i_file >> data)
		{
			Scale.push_back(stol(data)); // 第一个数是规模
			int all_num = Scale.back() * (Scale.back() + 1); //  当前规模行列式元素的总个数
			std::vector<float> temp_data;
			while(all_num--)
			{
				i_file >> data;
				temp_data.push_back(stof(data));
			}
			Data.push_back(temp_data);
		}
	}
	i_file.close();

	std::ofstream o_file;
	o_file.open(outputResult_path); // 清空计算结果文件
	o_file.close();
}

void saveOutputResult(buffer<float, 2> buf, int scale)
{
	std::ofstream o_file;
	o_file.open(outputResult_path, std::ios::app); //保存计算结果
	o_file << "问题规模 : " << scale <<std::endl;
	host_accessor m{ buf ,read_only };
	auto range = m.get_range();
	for (int i = 0; i < range[0]; i++) {
		for (int j = 0; j < range[1]; j++) {
			o_file << std::setw(16) << m[i][j];
		}
		o_file << std::endl;
	}
	o_file.close();
}

// 数据初始化(行列式)
void matrixInit(buffer<float, 2>& buf, std::vector<float> Data)
{
	host_accessor m{ buf ,read_write };

	auto range = m.get_range();
	int i_Data = 0;
	for(int i = 0; i < range[0]; i++)
	{
		for(int j = 0; j < range[1]; j++)
		{
			m[i][j] = Data[i_Data];
			i_Data += 1;
		}
	}
}

// 回代 把答案更新在最后一列
void backSubstitution(buffer<float, 2>& buf, queue& q)
{
	host_accessor m{ buf ,read_write }; //创建 主机访问器 将设备中数据拷贝到主机当中
	int n = m.get_range()[0];
	for (int i = n - 2; i >= 0; i--)
	{
		for (int j = n - 1; j - i >= 1; j--)
		{
			m[i][n] = m[i][n] - m[i][j] * m[j][n];
		}
	}
}

int getMax(buffer<float, 2>& buf, int k, int row_num)
{
	host_accessor m{ buf ,read_write }; //创建 主机访问器 将设备中数据拷贝到主机当中
	int max_index;
    float max;
	max_index = k;
    max= fabs(m[k][k]); // 把第一个要比较的数的绝对值作当前最大值
    for(int i = k + 1; i < row_num; i++)
    {
        if(max < fabs(m[i][k]))
        {
            max = fabs(m[i][k]);
            max_index = i;
        }
    }
	return max_index;
}

void colSwap(buffer<float, 2>& buf, int k, int max_index)
{
	host_accessor m{ buf ,read_write }; //创建 主机访问器 将设备中数据拷贝到主机当中
	int col_num = m.get_range()[1];
	// 若不是当前主元最大，则最大主元行与当前行交换
	for(int i = 0; i < col_num; i++)
	{
		float temp = m[k][i];
		m[k][i] = m[max_index][i];
		m[max_index][i] = temp;
	}
}

// 列主元消元法的并行算法
void GEColParallel(buffer<float, 2>& buf, queue& q)
{
	int n = buf.get_range()[0];
	for (int k = 0; k < n; k++)
	{
		int max_index = getMax(buf, k, n);

		if(max_index != k)
		{
			colSwap(buf, k, max_index);
		}

		q.submit([&](handler& h)
		{
			accessor m{ buf, h, read_write };  // 创建 缓冲区访问器 以访问设备上的缓冲区数据
			h.parallel_for(range(n - k) + 1, [=](auto idx)
			{
				int j = k + idx;
				m[k][j] = m[k][j] / m[k][k];
			}); // parallel_for卸载到加速器计算
		});

		q.submit([&](handler& h)
		{
			accessor m{ buf, h, read_write };
			h.parallel_for(range(n - (k + 1), n - k ), [=](auto idx)
				{
					int i = k + 1 + idx.get_id(0);
					int j = k + 1 + idx.get_id(1);
					m[i][j] += m[k][j] * (- m[i][k]);
				});
		});

		q.submit([&](handler& h)
		{
			accessor m{ buf, h, read_write };
			h.parallel_for(range(n - (k + 1)), [=](auto idx)
				{
					int i = k + 1 + idx;
					m[i][k] = 0;
				});
		});

	}
	q.wait();
}

// 求解过程
double test(int scale, std::vector<float> Data, const gauss_func& gauss, int times, queue& q)
{
	int row_num = scale;
	int col_num = scale + 1;
	buffer<float, 2> buf(range(row_num, col_num)); //创建一个scale * (scale + 1)大小的行列式 buffer 创建缓存区

	// 先运行一遍，以防止初始化操作、cache等其他因素影响算法执行时间的测量
	matrixInit(buf, Data); // 行列式初始化
	GEColParallel(buf, q); // 先运行一遍，以防止初始化操作、cache等其他因素影响算法执行时间的测量
	backSubstitution(buf, q);

	std::chrono::duration<double, std::milli> elapsed{}; // 时间记录变量
	for (int i = 0; i < times; i++)
	{
		// 行列式初始化
		matrixInit(buf, Data);
		auto start = std::chrono::high_resolution_clock::now();  // 运行时间记录
		GEColParallel(buf, q);
		backSubstitution(buf, q);
		auto end = std::chrono::high_resolution_clock::now();
		elapsed += end - start;
	}
	saveOutputResult(buf, scale);
	return elapsed.count() / times;
}

// 输出结果(时间)
void testOutput(const std::vector<gauss_func>& gauss_funcs,
	std::string name,
	int times,
	queue& q)
{

	std::vector<int> Scale;
	std::vector<std::vector<float>> Data; // 每一个Data[i]里面存着Scale[i]规模的一个顺序展开的行列式
	dataInit(Scale, Data); // 把待计算的行列式加载到内存

	std::cout << std::endl;
	std::cout << "  ******* " << name << " *******  " << std::endl;
	std::cout << std::endl;
	std::cout << "问题规模    耗时";
	std::cout << std::endl;

	for (int i = 0; i < Scale.size(); i++)
	{
		std::cout << Scale[i] << "           ";
		for (auto& func : gauss_funcs)
		{
			std::cout << test(Scale[i], Data[i], func, times, q) << " ";
		}
		std::cout << std::endl;
	}
}

int main() 
{

	gpu_selector selector;
	queue q(selector);
	device my_device = q.get_device();

	std::cout << "Device: " << my_device.get_info<info::device::name>() << std::endl;

	std::vector<gauss_func> gauss_funcs =
	{
		GEColParallel,
	};

	std::string name = "列主元消元法---并行";
	const int times = 10;

	testOutput(gauss_funcs, name, times, q);
}
